Method and apparatus for transforming uniformly or non-uniformly sampled interferograms to produce spectral data

ABSTRACT

A reconstruction matrix used for calculating a hyperspectral data-cube includes rows of periodic functions. Each row of the reconstruction matrix corresponds to a selected wavelength and each column corresponds to a selected retardance of an interferometer. The periodic functions have as a parameter the selected wavelength of the corresponding row and are sampled at the selected retardances of each of the corresponding columns. An interferogram data-cube is obtained and includes an array of one or more simultaneously measured interferograms. Each row of the interferogram data-cube corresponds to one of the selected retardances and each column corresponds to a different interferogram from the simultaneously measured interferograms. A set of matrix-vector products for each of the interferograms is formed by multiplying the reconstruction matrix with a column of the interferogram data-cube to form the hyperspectral data-cube.

SUMMARY

The present disclosure is directed to a method and apparatus fortransforming uniformly or non-uniformly sampled interferograms toproduce spectral data at predefined wavelengths. In one embodiment, aset of reference retardances of an interferometer are determined, as area set of wavelengths corresponding to spectral slices of a hyperspectraldata-cube. A reconstruction matrix is formed that includes rows ofperiodic functions. Each row of the reconstruction matrix corresponds toa selected wavelength of the set of wavelengths, and each column of thereconstruction matrix corresponds to a selected retardance of thereference retardances. The periodic functions have as a parameter theselected wavelength of the corresponding row and are sampled at theselected retardances of each of the corresponding columns. Aninterferogram data-cube is obtained that includes an array of one ormore simultaneously measured interferograms. Each row of theinterferogram data-cube corresponds to one of the selected retardancesand each column of the interferogram data-cube corresponds to adifferent interferogram from the simultaneously measured interferograms.A set of matrix-vector products is formed for each of theinterferograms. Each of the matrix-vector products includes a matrixmultiplication of the reconstruction matrix with a column of theinterferogram data-cube. The set of matrix-vector products forms thehyperspectral data-cube.

These and other features and aspects of various embodiments may beunderstood in view of the following detailed discussion and accompanyingdrawings.

BRIEF DESCRIPTION OF THE DRAWINGS

The discussion below makes reference to the following figures, whereinthe same reference number may be used to identify the similar/samecomponent in multiple figures. The drawings are not necessarily toscale.

FIG. 1 is a diagram of an apparatus according to an example embodiment;

FIG. 2 is a set of plots showing interferograms processed via an opticaldevice according to example embodiments;

FIG. 3 is a graph of retardance versus time for an optical deviceaccording to an example embodiment;

FIG. 4 is a set of plots showing the measurement of a wavelengthdependence of a zero path delay of an interferometer according to anexample embodiment;

FIG. 5 is a diagram of a reconstruction matrix according to an exampleembodiment;

FIG. 6 is a conceptual diagram of an interferogram data-cube accordingto an example embodiment;

FIG. 7 is a diagram of a matrix that represents an interferogramdata-cube according to an example embodiment;

FIG. 8 is a diagram of a system response matrix according to an exampleembodiment;

FIG. 9 is a set of plots showing the effect of ramp and window functionson a reconstruction matrix according to example embodiments; and

FIG. 10 is a flowchart of a method according to an example embodiment.

DETAILED DESCRIPTION

The present disclosure relates to signal processing in hyperspectralimaging systems. Hyperspectral imaging systems described herein use apolarization interferometer that is configured to introduce a variableoptical path delay (or retardance) in components of light that passthrough the interferometer. The device that causes the path delay(referred to as a variable optical retarder) is placed between twopolarizers such that the variable path delay is introduced between firstrays in an incident polarization direction and second rays in anorthogonal polarization (e.g., ordinary and extraordinary rays), withboth sets of rays tracing a common path. This path delay causes awavelength-dependent phase shift between the first and second rays. Thepath delay causes light exiting the polarization interferometer to forminterferograms that are detected via an optical sensor, e.g., afocal-plane array.

A polarization interferometer can use one or more liquid-crystal (LC)cells as a variable optical retarder. Such a variable optical retarderis referred to herein as a liquid-crystal variable retarder (LCVR).Generally, liquid-crystal (LC) materials are liquids having somecrystalline properties (e.g., orientation of internal structures, suchas the LC director that indicates the local average alignment of LCmolecules) that can be selectably altered by applying an externalstimulus, such as an electric field or a magnetic field. A change inorientation of the LC director alters the optical properties of the LCmaterials, e.g., changing the optical axis of the LC birefringence.

An LCVR generates a variable optical path delay, or a variableretardance, between two orthogonal polarizations of light that travelthrough the liquid crystal. One or more liquid-crystal cells within theLCVR function as electrically tunable birefringent elements. By varyingthe voltage across the electrodes of the liquid-crystal cell, the cellmolecules change their orientation, and it is possible to controllablychange the optical path delay over a period of time.

To create a polarization interferometer with an LCVR, the LCVR is placedbetween a first polarizer and a second polarizer with nominally parallelor perpendicular polarization axes. The slow axis of the LCVR (thepolarization axis with the variable optical path delay) is orientednominally 45 degrees with respect to the polarization direction of thefirst polarizer. Incoming light is polarized to an incident polarizationdirection by the first polarizer. Because the slow axis of the LCVR isat 45 degrees with respect to this incident polarization direction, thepolarized incident light can be described in terms of a portion of lightpolarized parallel to the slow axis of the LCVR and a portion of lightpolarized perpendicular to this axis.

As the light passes through the LCVR, it acquires a wavelength-dependentrelative phase shift between the first and second polarizations, therebyleading to a wavelength-dependent change in the polarization state. Thesecond polarizer, or analyzer, oriented either parallel or perpendicularto the first polarizer, interferes the portion of light polarizedparallel to the slow axis of the LCVR with the portion of lightpolarized perpendicular, changing the wavelength-dependent polarizationstate at the output of the LCVR into a wavelength-dependent intensitypattern that can be sensed by an optical detector or a focal planearray. By sensing this intensity while varying the retardance of theLCVR, it is possible to measure an interferogram of the incoming light,which can be used to ascertain spectral properties of the incominglight.

As noted above, a polarization interferometer is used in hyperspectralimaging applications because of its abilities to encode spectralinformation of the incident light into an intensity pattern that iseasily measured with a non-spectrally-resolving detector. Hyperspectralimaging refers to methods and devices for acquiring hyperspectraldatasets or data-cubes, which may include images where densely sampled,finely resolved spectral information is provided at each pixel.

The wavelength-dependent intensity pattern provided by the polarizationinterferometer corresponds approximately to a cosine transform of thespectrum of the incident light. By recording each pixel's intensitypattern at the output of a polarization interferometer as a function ofthe LCVR's retardance, the interferograms generated by all points of ascene imaged through the LCVR can be sampled simultaneously. From this,the hyperspectral data-cube can be nominally recovered by applying atransform, such as an inverse cosine transform or Fourier transformalong the retardance axis, to the recorded spatially-dependentinterferogram.

Hyperspectral imagers based on liquid-crystal polarizationinterferometers generate interferograms that may be non-uniformlysampled with respect to optical path delay. A standard Fourier transformmay not be optimal for converting the raw data generated by such adevice to a hyperspectral image. This disclosure describes methods andapparatuses used to reconstruct the hyperspectral data-cube with betteraccuracy and fewer artifacts than with a Fourier transform, whether ornot the interferograms are uniformly sampled.

To enhance understanding of hyperspectral image processing describedherein, a block diagram in FIG. 1 illustrates an apparatus 100 thatperforms image processing according to an example embodiment. Theapparatus 100 includes a device controller 102, which may include one ormore processors, such as central processing units, subprocessors,graphics processing units, digital signal processors, etc. Thecontroller 102 is coupled to a memory 104 that includes functionalmodules that will be described in greater detail below. The memory 104may include a combination of volatile and non-volatile memory, and maystore instructions and data as known in the art.

The apparatus 100 includes an optical section 106 with an externaloptical interface 108 that receives light from outside the apparatus100. The external optical interface 108 may include windows, lenses,filters, apertures, etc., suitable for passing light 109 from outsidethe apparatus 100 to internal optical components. In this example, theexternal optical interface 108 is shown coupled to an external lens 110.

A polarization interferometer 112 is located in the optical section 106of the apparatus 100. The polarization interferometer 112 is coupled tothe controller 102, e.g., via electrical signal lines. The controller102 applies signals to the polarization interferometer 112 to cause atime-varying optical path delay or retardance in an LCVR 112 a that ispart of the interferometer 112. This time-varying optical path delaycauses a wavelength-dependent phase shift between differentpolarizations of the light 109, resulting in interferograms that vary asa function of the optical path delay. The interferograms are detected byan image sensor 114 (e.g., an array of sensor pixels, focal plane array)which is also coupled to the controller 102.

A retardance controller 118 instructs the device controller 102 to applya control signal to the LCVR 112 a to achieve a time-varying retardancetrajectory. An image processor 120 uses this retardance trajectory as ameasure of time-varying path delay together with interferograms detectedat the image sensor 114. Each detected interferogram can be processed bycalculating a transform as a function of the path delay at acorresponding position of the LCVR 112 a, and together the processedinterferograms comprise spatially-dependent spectral data, e.g., ahyperspectral data-cube.

All the processing described herein can take place within the apparatus100 or via external computer 124 coupled by interface 122 or somecombination thereof. This is indicated by image processor 126 operatingon computer 124. The computer 124 and/or the apparatus 100 may include agraphics processing unit (GPU) 128, 130 that can be used to speed upcertain computations as described below. Generally, a GPU 128, 130 willhave a large number of processing cores (e.g., hundreds or thousands)configured for parallelizable sub-computations, such as matrix orsub-matrix multiplications. Assuming each of the sub-computations areindependent, a GPU can greatly reduce image processing time. As thecalculations below involve parallelizable matrix multiplication, a GPUmay be used to improve processing speeds.

Hyperspectral imagers described herein generate raw signals that areinterferograms with respect to optical path delay independently at eachpixel. When these systems are based on liquid-crystal polarizationinterferometers, the sampling of the interferograms may not be entirelyuniform with respect to optical path delay. As a result, typical FastFourier Transform (FFT) algorithms may not work to calculatehyperspectral data-cubes from the interferograms. Furthermore, even ifthe interferogram sampling were uniform, an FFT of the interferogramsmay not lead to an optimal estimate of spectral data, e.g., thehyperspectral data-cube. Non-uniform sampling may be intentionallyintroduced, for example through denser sampling at shorter path delays,in order to enhance the sensitivity at certain wavelengths.

One way that non-uniform sampling might be remedied prior to applying anFFT to an interferogram is to resample the interferogram at uniformintervals of optical path delay or retardance via interpolation. One canassume that the optical path delay of each recorded interferogram sampleis known, for example by obtaining the phase of an interferogram of amonochromatic reference light source of known wavelength as a functionof time. If the optical path delay were changed monotonically as afunction of time as is generally the case, there would be a one-to-onerelationship between optical path delay and time. Individualinterferogram samples would be taken at known retardances, and auniformly spaced retardance grid could be used to calculate, viainterpolation of the interferogram samples, the interpolatedinterferograms.

However, this interpolation method may be computationally expensive,especially if higher-order (e.g., cubic spline) methods are used orimage sizes are large. Also, in regions of the interferogram that may besampled near the Nyquist limit or that may have few samples perinterference fringe (especially for shorter wavelengths near zero pathdelay), the interpolation may cause significant loss and/or distortionof the signal, which would be disadvantageous and potentially causeartifacts in the reconstructed spectral data. The hyperspectral imagersdescribed herein may already be less sensitive to shorter wavelengths,and such an interpolation just exacerbates this problem.

In FIG. 2, a set of plots show how different sampling methods can affectthe signal strength and distortion of an interferogram according to anexample embodiment. The top plot 200 is an interferogram of a blue LEDmeasured through a polarization interferometer based on an LCVR. Theinterferogram was sampled at uniformly spaced time intervals while theretardance of the LCVR was scanned at a nominally constant rate.However, due to small fluctuations in the scanning rate of the LCVR, theinterferogram in the top plot 200 was sampled at non-uniform intervalsof retardance. The middle plot 202 shows the results of resampling, viainterpolation, of the interferogram in the top plot 200 at uniformretardance increments. This interpolation decreases the interferogramsignal strength (as measured by the interferogram's standard deviation)by almost 20%. Plot 204 shows a separately measured interferogram of thesame blue LED that was purposely sampled non-uniformly with respect toretardance (although at uniformly-spaced time intervals, as shown inFIG. 3) to increase the number of samples near zero path delays 206-208which is where, as can be seen in plots 200, 202, 204, most of thesignal from short (blue) wavelengths is located. In this example, theintentional non-uniform sampling of the interferogram in 204 boosted thestandard deviation of the signal with respect to the interferogram in200 by about 30%. In general, contributions from wavelengths across thedetectable optical spectrum (e.g., 400 nm-1000 nm for CMOS detectors)will be present in measured interferograms. The intentional non-uniformretardance sampling described above can have the effect of increasingthe signal strength at shorter wavelengths with respect to longerwavelengths, which is often desirable to maintain dynamic range acrossthe detectable optical spectrum.

The non-uniform retardance sampling can also be dealt with via a numberof methods developed for Fourier transforms of non-uniformly sampleddata, such as the non-uniform discrete Fourier transform, which can becalculated with the NuFFT software library. However, a Fourier transformacross all frequencies may be overkill, as many of the calculatedfrequencies may not even correspond to wavelengths within the detectableoptical spectrum. Instead, the image processor 120 can approximate aFourier transform by taking a dot-product of the interferogram with onlythe (uniformly or non-uniformly sampled) Fourier components of interest,or equivalently, only those rows of the (regular or non-uniform)discrete Fourier transform (DFT) matrix that contribute to wavelengthsof interest.

An advantage of this method is that each wavelength can be treatedindependently to incorporate different apodizations andwavelength-dependent zero path delay, and the Fourier components can besampled with the same non-uniformity as the interferogram such that thenon-uniformity is no longer an issue. Any additional lineartransformations can also be combined with this matrix to form an overall“reconstruction” matrix for directly reconstructing spectral data fromthe raw interferograms. Such a matrix could even have a spatial index toincorporate spatial calibrations such as a spatially-dependentretardance.

To generate the reconstruction matrix, a reference retardance vs. timecurve is obtained from a sensor, as shown in by the example curve 300 inFIG. 3. For example, a monochromatic light source 116 (see FIG. 1)located on one side of a liquid-crystal cell may be used to obtain thiscurve 300. The light from the light source 116 is detected by the sensor114 while the controller 102 applies a voltage to the LC cell(s) of theLCVR 112 a, causing the retardance of the LCVR to transition from afirst to a second retardance. The controller 102 then analyzes thedetected oscillations in intensity of the light source 116 to determinethe time-dependent retardance, or equivalently, the phase of theinterferogram of the light source 116 for each sample of theinterferogram. Other ways of obtaining the data of plot 300 compriseanalyzing a capacitance measurement of the LC cell(s) using a similarprocedure.

The system may have wavelength- and temperature-dependent variations inzero path delay which may arise, for example, by compensating a residualliquid-crystal retardance with another birefringent material. It isdesirable that interferograms are measured at both positive and negativepath delays, and in some cases, the residual liquid-crystal retardanceprevents this and must be compensated. One way to compensate theresidual retardance is to include within the polarization interferometera fixed retardance of opposite sign, e.g., a multiple-order waveplatewith slow axis perpendicular to the slow axis of the LCVR. If thiscompensating waveplate has a retardance that depends on wavelengthand/or temperature in a manner different than the retardance of theLCVR, the zero path delay of the polarization interferometer woulddepend on both wavelength and temperature. In such a case, the systemcan also determine and store a zero path delay point, with respect tothe reference retardance curve, for each wavelength and thus each row ofthe reconstruction matrix. An example of determining zero path delayaccording to an example embodiment is shown in a set of plots in FIG. 4.

In FIG. 4, plot 400 is a set of measurements of interferograms of anarrowband light source stepped in 10 nm increments, each shown stackedupon the other. Wavelength increases up the vertical axis. Plot 410 isan expanded view of the portion 402 of plot 400, and plot 420 is asubset of the portion 412 of plot 410. The highlighted points 422 inplot 420 indicate the reference retardance value that represents theactual zero path delay for each wavelength, σ(λ)=0, where σ is theactual interferometer path delay and not the measured referenceretardance. These points 422 form a curve and not a vertical linebecause some of the liquid-crystal retardance is compensated with awaveplate of different dispersion, thus resulting in a residualdispersion across the detectable optical spectrum.

This residual dispersion is related to a sum of indices of refraction ofboth the liquid-crystal material and the compensation waveplate. Asindicated by plot 430 and equation 432, a Cauchy approximation (e.g., upto 6th order) can be used to fit the data. Note that the coefficients ofthe equation generally depend on temperature. Using equation 432, thereference retardance corresponding to zero path delay can be calculatedat any wavelength where the Cauchy approximation is valid. This equationcan also be solved for different temperatures, assuming the coefficientvariation with respect to temperature has been characterized. Thecoefficients can be stored in memory of the apparatus 100 as known inthe art and used by the image processor in the computations describedbelow.

Another parameter that can be used by the image processor 120 is therate of change of phase delay of a given wavelength of interest (e.g.,λ_(i), corresponding to the i^(th) row of the reconstruction matrix)relative to the rate of change of phase of the reference retardancecurve at a reference wavelength λ₀. A scaling factor can be obtained,for example, by knowledge of the liquid crystal indices of refractionn_(e), n₀ as a function of wavelength λ and temperature T, where the e-and o-subscripts represent the respective indices for the extraordinaryand ordinary rays. Equation (1) below can be used to determine thescaling factor α_(i), at wavelength λ_(i), or a calibration can be usedto empirically determine a value.

$\begin{matrix}{\alpha_{i} = {\frac{{n_{e}\left( {\lambda_{i},T} \right)} - {n_{o}\left( {\lambda_{i},T} \right)}}{{n_{e}\left( {\lambda_{0},T} \right)} - {n_{o}\left( {\lambda_{0},T} \right)}}\left( \frac{\lambda_{0}}{\lambda_{i}} \right)}} & (1)\end{matrix}$

In view of the above, the elements A_(ij) of the reconstruction matrixmay then be represented as shown in Equation (2) below, where Φ_(j) isthe phase of the reference retardance curve corresponding tointerferogram sample j. The phase Φ_(j) is found empirically byanalyzing the interferogram of a reference light source, e.g., usingTakeda's method; it relates to the reference retardance Γ_(j) asΦ_(j)=2πΓ^(j)/λ₀ The angle ϕ_(i) is the phase of the referenceretardance curve corresponding to zero path delay for wavelength λ_(i),calculated e.g. as in Equation (3) via the fitting coefficients obtainedin FIG. 4. Subtracting the phase ϕ_(i) from the phase of the referenceretardance curve Φ_(j) phase shifts the periodic function sampled in rowi such that it has a phase of 0 at an interferogram sample correspondingto the zero-retardance point of the interferometer at wavelength λ_(i).

$\begin{matrix}{{{A_{ij} = {\cos\left( {\alpha_{i}\left( {\Phi_{j} - \phi_{i}} \right)} \right)}};\mspace{25mu}{i \in \left\{ {1\mspace{14mu}.\;.\;.\mspace{14mu} M} \right\}}},\mspace{14mu}{j \in \left\{ {1\mspace{14mu}.\;.\;.\mspace{14mu} N} \right\}}} & (2) \\{\phi_{i} = {2\;{{\pi\left( {B + \frac{C}{\lambda_{i}^{2}} + \frac{D}{\lambda_{i}^{4}} + \frac{E}{\lambda_{i}^{6}}} \right)}/\lambda_{0}}}} & (3)\end{matrix}$

Note the reconstruction matrix elements in Equation (2) only need to becalculated once per image acquisition, as each acquisition has aspecific reference retardance curve and temperature. In FIG. 5, adiagram shows a reconstruction matrix 500 according to an exampleembodiment. The matrix 500 includes rows of periodic functions 502(e.g., cosine functions) as in Equation (2). Each row corresponds to aselected wavelength of a set of wavelengths 504, which each correspondto spectral slices of a hyperspectral data-cube. The wavelengths 504 canbe preselected to have a particular range and distribution. Note thatthe previously described reference wavelength λ₀ may or may not beincluded in the set of wavelengths 504.

Each column of the matrix 500 corresponds to a selected retardance of aset of reference retardances 506. The set of reference retardances 506are associated with a state of an interferometer at a particular pointin time. Note that the set of retardances 506 may be distributed atnon-uniform retardance intervals (e.g., as shown in FIG. 3), and theinitial retardance Γ₁ may or may not be a minimal retardance.

The periodic functions 502 have as a parameter the selected wavelengths504 (or the indices of the selected wavelengths 504) of thecorresponding row and are sampled at the reference retardances 506 ofeach of the corresponding columns. The matrix 500 can be multiplied withan interferogram data-cube to obtain spectral data. In FIG. 6, aconceptual diagram shows an interferogram data-cube according to anexample embodiment. Elements 600 represent individual sensors of adetector, e.g., pixels of a focal plane array, arranged across anxy-plane. The t-axis represents a change in retardance of a polarizationinterferometer from Γ₁ to Γ_(N), which each correspond to a differentpoint in time between t₁ and t_(N) during which the data-cube isacquired. The curves 602 represent individual interferograms measured ateach sensor 600 over the acquisition time period. Therefore, theinterferogram data-cube is formed of an array with each element having asample value corresponding to a particular sensor position andretardance value. Points 604 represent individual samples correspondingto retardance Γ_(j) for the sensors 600.

While the illustrated sensors 600 are in a two-dimensional spatialarrangement, for purposes of calculation they can be labeled by a singleindex, e.g., p_(k), k ∈{1 . . . P}. As such, the interferogram data-cubecan be represented as a two-dimensional matrix X 700 as shown in FIG. 7.The matrix 700 includes samples 702 that correspond to simultaneouslymeasured interferograms. Each row of the matrix 700 corresponds to aselected retardance of the reference retardances 506 and each columncorresponds to a different interferogram measurement from the set ofsimultaneously measured interferograms with index k. As shown here, theinterferogram samples are labeled by sensor positions 704, the sensorarray having P positions, e.g., P=number of pixel rows* number of pixelcolumns.

The calculation of spectral data or the hyperspectral data-cube involvesforming a set of matrix-vector products for each of the interferograms,each of the matrix-vector products comprising a matrix multiplication ofthe matrix 500 with a column k of the interferogram data-cube 700. Thus,the resulting hyperspectral data-cube H can be formed by the matrixmultiplication H=AX. Each row of H can be composed into an imagecorresponding to one of the wavelengths 504, each pixel of the imagecorresponding to an intensity at the corresponding wavelength and at thespatial location labeled with index k.

In some embodiments, different reconstruction matrices 500 can beconstructed for at least two spatial regions of the interferometer.Thus, different pixel locations p_(k) in FIG. 7 can be multiplied by adifferent reconstruction matrix. This is indicated in FIG. 7 by regions705-707, which may each be associated with a different reconstructionmatrix 500. Such regions may encompass one or more individual pixels andmay also encompass non-contiguous portions of the data-cube matrix 700.In this way, compensation for spatially-dependent characteristics of thepolarization interferometer (or other optical components) can be made atany desired granularity.

The reconstruction of spectral data from interferograms that areasymmetric with respect to zero path delay requires special handling toprevent the formation of low-resolution spectral artifacts (see, e.g.,C. D. Porter and D. B. Tanner, “Correction of phase errors in Fourierspectroscopy,” Int. J. Infrared Millimeter Waves 4, 273-298 (1983)). Oneway to incorporate this special handling into the reconstruction matrix500 is to multiply each row of the matrix 500 with a linear “ramp”function that starts equal to zero at the end nearest zero path delay,goes through one half at the (wavelength-dependent) zero path delaypoint, and continues to one with the same slope, at which point itremains constant. Without loss of generality, it can be assumed that thephase of the reference retardance curve starts at a minimum of 0 forΦ_(j,j=1) and increases monotonically for increasing sample number j.Then, the wavelength-dependent ramp can be defined as in Equation (4)below. Equation (4) defines a matrix that can be multiplied elementwisewith the reconstruction matrix to accommodate asymmetric interferograms.

$\begin{matrix}{r_{ij} = \left\{ {{\begin{matrix}{\frac{\Phi_{j}}{2\;\phi_{i}},} & {\Phi_{j} < {2\phi_{i}}} \\{1,} & {\Phi_{j} \geq {2\phi_{i}}}\end{matrix};\mspace{20mu}{i \in \left\{ {1\mspace{14mu}.\;.\;.\mspace{14mu} M} \right\}}},\mspace{14mu}{j \in \left\{ {1\mspace{14mu}.\;.\;.\mspace{14mu} N} \right\}}} \right.} & (4)\end{matrix}$

The rows of the reconstruction matrix 500 can also be multiplied by awindowing or apodization function as long as it is symmetric about zeropath delay for each wavelength. One advantage of calculating thereconstruction matrix row-by-row is that a different size windowingfunction can be applied for each wavelength. For example, due tointrinsic path delay inhomogeneities across an LCVR that might limit thebest obtainable spectral resolution for shorter wavelengths, shorterwavelengths typically generate interferograms with most of their signalenergy located around zero path delay. Therefore, when reconstructingthe shorter wavelength portions of a spectrum from an interferogram, itis advantageous to first multiply the interferogram by a windowingfunction that includes just those portions where there is signal fromshorter wavelengths. Making the window bigger would add more noise butnot signal, as unlike the signal, the noise is often distributed acrossall path delays. And, to avoid repeating multiplications of theinterferograms at different locations with the same window function(s),the window functions can be multiplied elementwise with the rows of thereconstruction matrix 500.

If A is the reconstruction matrix 500 and B is a matrix with columnsthat are the system response at each wavelength represented by a row ofthe reconstruction matrix 500, then the reconstruction matrix 500 can bepre-multiplied by the diagonalization matrix (AB)⁻¹ in order to“diagonalize” the reconstruction. That is, the full reconstructionmatrix could be set to (AB)⁻¹ A so that multiplication by the systemresponse matrix B (or by the measured response at a given wavelength)produces the identity. This is analogous to a least-squares solution tothe reconstruction of spectral data from interferograms. Also, becauseAB may be ill-conditioned, it may be necessary to provide some form ofregularization to prevent the above matrix inversion from diverging. Forexample, one could add a multiple of the identity prior to inversion:(AB+βI)⁻¹, where I is the M×M identity matrix and β is a multiplicationfactor. The reconstruction matrix 500 can also include information aboutthe spectral sensitivity of the interferometer in order to equalize thesystem response across wavelengths. The diagram in FIG. 8 shows thesystem response matrix B 800. The columns of the matrix 800 are theexpected raw interferograms generated by a unit wavelength stimulus tothe interferometer at the wavelengths 504 corresponding to each of therows of A.

In FIG. 9, a pair of stacked plots shows two reconstruction matrices asdescribed above and plotted similarly to the measurements in FIG. 4.Each line of plot 900 is a row of the reconstruction matrix A as inEquation (2), accounting for non-uniform retardance spacing,wavelength-dependent zero path delay, and liquid-crystal dispersion. Theplot 902 shows the effects of additional processing, which includesmultiplication of each row by the ramp function of Equation (4) to avoiddouble-counting the symmetric part of the interferogram, as well aswavelength-dependent windowing to separately define the resolution binsize for each wavelength.

In FIG. 10, a flowchart shows a method according to an exampleembodiment. The method involves determining 1000 a set of referenceretardances of an interferometer. A set of wavelengths corresponding tospectral slices of a hyperspectral data-cube are also determined 1001,e.g., based on a desired spectral range of the hyperspectral data-cube.A reconstruction matrix is formed 1002. The reconstruction matrix hasrows of sampled periodic functions. Each row of the reconstructionmatrix corresponds to a selected wavelength of the set of wavelengths,and each column of the reconstruction matrix corresponds to a selectedretardance of the reference retardances. The periodic functions have asa parameter the selected wavelength of the corresponding row and aresampled at the reference retardances of each of the correspondingcolumns.

An interferogram data-cube is obtained 1003, e.g., via an opticaldetector. The interferogram data-cube includes an array of one or moresimultaneously measured interferograms. The simultaneously measuredinterferograms may be non-uniformly sampled with respect to optical pathdelay. Each row of the interferogram data-cube corresponds to one of theselected retardances and each column of the interferogram data-cubecorresponds to a different interferogram from the set of simultaneouslymeasured interferograms. A set of matrix-vector products is formed 1004for each of the interferograms. Each of the matrix-vector productsincludes a matrix multiplication of the reconstruction matrix with acolumn of the interferogram data-cube. The set of matrix-vector productsforms the hyperspectral data-cube.

The various embodiments described above may be implemented usingcircuitry, firmware, and/or software modules that interact to provideparticular results. One of skill in the relevant arts can readilyimplement such described functionality, either at a modular level or asa whole, using knowledge generally known in the art. For example, theflowcharts and control diagrams illustrated herein may be used to createcomputer-readable instructions/code for execution by a processor. Suchinstructions may be stored on a non-transitory computer-readable mediumand transferred to the processor for execution as is known in the art.The structures and procedures shown above are only a representativeexample of embodiments that can be used to provide the functionsdescribed hereinabove.

Unless otherwise indicated, all numbers expressing feature sizes,amounts, and physical properties used in the specification and claimsare to be understood as being modified in all instances by the term“about.” Accordingly, unless indicated to the contrary, the numericalparameters set forth in the foregoing specification and attached claimsare approximations that can vary depending upon the desired propertiessought to be obtained by those skilled in the art utilizing theteachings disclosed herein. The use of numerical ranges by endpointsincludes all numbers within that range (e.g. 1 to 5 includes 1, 1.5, 2,2.75, 3, 3.80, 4, and 5) and any range within that range.

The foregoing description of the example embodiments has been presentedfor the purposes of illustration and description. It is not intended tobe exhaustive or to limit the embodiments to the precise form disclosed.Many modifications and variations are possible in light of the aboveteaching. Any or all features of the disclosed embodiments can beapplied individually or in any combination are not meant to be limiting,but purely illustrative. It is intended that the scope of the inventionbe limited not with this detailed description, but rather determined bythe claims appended hereto.

What is claimed is:
 1. A method, comprising: determining a set ofreference retardances of an interferometer; determining a set ofwavelengths corresponding to spectral slices of a hyperspectraldata-cube; forming a reconstruction matrix comprising rows of periodicfunctions, wherein each row of the reconstruction matrix corresponds toa selected wavelength of the set of wavelengths, and each column of thereconstruction matrix corresponds to a selected retardance of thereference retardances, wherein the periodic functions have as aparameter the selected wavelength of the corresponding row and aresampled at the selected retardance of each of the corresponding columns;obtaining an interferogram data-cube comprising an array of one or moresimultaneously measured interferograms, wherein each row of theinterferogram data-cube corresponds to one of the selected retardancesand each column of the interferogram data-cube corresponds to adifferent interferogram from the simultaneously measured interferograms;and forming a set of matrix-vector products for each of theinterferograms, each of the matrix-vector products comprising a matrixmultiplication of the reconstruction matrix with a column of theinterferogram data-cube, the set of matrix-vector products forming thehyperspectral data-cube.
 2. The method of claim 1, wherein theinterferograms correspond to different image points within theinterferometer.
 3. The method of claim 1, wherein the simultaneouslymeasured interferograms are non-uniformly sampled with respect tooptical path delay.
 4. The method of claim 1, wherein the interferometercomprises a polarization interferometer.
 5. The method of claim 4,wherein the polarization interferometer comprises a liquid-crystalvariable retarder, and wherein a relative rate of oscillation of theperiodic functions in each row is based on a liquid-crystalbirefringence of the liquid-crystal variable retarder as a function ofat least one of wavelength and temperature.
 6. The method of claim 1wherein the calculation takes place on a graphics-processing unit. 7.The method of claim 1 wherein different matrices are constructed for atleast two spatial regions of the interferometer.
 8. The method of claim1 wherein the periodic function is a sinusoidal function.
 9. The methodof claim 1 wherein each of the periodic functions is phase-shifted foreach of the selected wavelengths so that a phase of 0 correspondsroughly to a zero-retardance point of the interferometer at each of theselected wavelengths.
 10. The method of claim 9 wherein the phase shiftdepends on a temperature of the interferometer.
 11. The method of claim1, further comprising multiplying each row of the reconstruction matrixby a wavelength-dependent function to reduce artifacts caused byasymmetrically measured interferograms.
 12. The method of claim 1,wherein each row of the reconstruction matrix is multiplied by awavelength-dependent windowing or apodization function.
 13. The methodof claim 1, wherein each row of the reconstruction matrix is multipliedby the reciprocal of a sensitivity function at each selected wavelength.14. The method of claim 1, wherein the reconstruction matrix ismultiplied by a diagonalization matrix (AB)⁻¹ to yield (AB)⁻¹A, where Ais the reconstruction matrix and B is a system response matrix withcolumns that are the expected raw interferograms generated by a unitwavelength stimulus to the interferometer at the wavelengthscorresponding to each of the rows of A.
 15. The method of claim 14,wherein a regularization technique is applied to prevent thediagonalization matrix from diverging.
 16. An apparatus, comprising: aninterferometer operable to transition through a set of referenceretardances; an optical sensor operable to simultaneously measureinterferograms corresponding to different positions of the opticalsensor as the interferometer transitions through the set of thereference retardances; and a controller coupled to the interferometerand the optical sensor, the controller operable to: determine a set ofwavelengths corresponding to spectral slices of a hyperspectraldata-cube; form a reconstruction matrix comprising rows of periodicfunctions, wherein each row corresponds to a selected wavelength of theset of wavelengths, and each column corresponds to a selected retardanceof the reference retardances, wherein the periodic functions each haveas a parameter the selected wavelength of the corresponding row and aresampled at the selected retardance of the corresponding column; obtainan interferogram data-cube comprising an array of one or more of thesimultaneously measured interferograms, wherein each row of theinterferogram data-cube corresponds to one of the selected retardancesand each column of the interferogram data-cube corresponds to adifferent interferogram from the simultaneously measured interferograms;and form a set of matrix-vector products for each of the interferograms,each of the matrix-vector products comprising a matrix multiplication ofthe reconstruction matrix with a column of the interferogram data-cube,the set of matrix-vector products forming the hyperspectral data-cube.17. The apparatus of claim 16, wherein each of the periodic functions isphase-shifted for each of the selected wavelengths so that a phase of 0corresponds roughly to a zero-retardance point of the interferometer ateach of the selected wavelengths.
 18. The apparatus of claim 16, whereinthe controller is further configured to multiply each row of thereconstruction matrix by a wavelength-dependent function to reduceartifacts caused by asymmetrically measured interferograms.
 19. Theapparatus of claim 16, wherein each row of the reconstruction matrix ismultiplied by a wavelength-dependent windowing or apodization function.20. The apparatus of claim 16, wherein the reconstruction matrix ismultiplied by a diagonalization matrix (AB)⁻¹ to yield (AB)⁻¹A, where Ais the reconstruction matrix and B is a system response matrix withcolumns that are the expected raw interferograms generated by a unitwavelength stimulus to the interferometer at the wavelengthscorresponding to each of the rows of A.